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1 Introduction 



Neutron stars are among the most fascinating bodies in our universe. They 
contain over a solar mass of matter within a radius of ~ 10 km at densities of order 
10^^ g/cc. They probe the properties of cold matter at extremely high densities, 
and have proven to be fantastic test bodies for theories of general relativity. 
In a broader perspective, neutron stars and heavy ion collisions provide access 
to the phase diagram of matter at extreme densities and temperatures, that is 
basic for understanding the very early Universe and several other astrophysical 
phenomena. 

The discovery of the neutron by Chadwick in 1932 prompted Landau Q to 
predict the existence of neutron stars. The first theoretical calculations of neutron 
stars were performed by Oppenheimer and Volkoff in 1939 assuming that they 
are gravitationally bound states of neutron Fermi gas. The calculated stars had 
a maximum mass of ~ 0.7 M©, central densities up to ~ 6 x 10^^ g/cm^ and 
radii ~ 10 km. For comparison the density of nuclear matter inside a large 
nucleus like ^osp]-, jg ^ Qjg nucleons/fm^, i.e. ~ 2.7 xlO^^ g/cm^ Their 
predicted maximum mass was less than the Chandrasekhar mass limit of ~ 1.4 
Mq for white dwarfs made up of iron group nuclei, and having densities up to 
~ 10^ g/cm^ [^. The pressure to balance the gravitational attraction in white 
dwarfs and Oppenheimer- Volkoff neutron stars is supplied by degenerate electron 
and neutron Fermi gases respectively. 

In 1934 Baade and Zwicky |^] suggested that neutron stars may be formed in 
supernovae in which the iron core of a massive star exceeds the Chandrasekhar 
limit and collapses. The large amount of energy released in the collapse blows 
away the rest of the star and the collapsed core may form a neutron star. For 
efficient production of neutron stars with this mechanism, the maximum mass 
of neutron stars should exceed 1.4 Mq. In the 60's, using schematic models of 
nuclear forces, Tsuruta and Cameron showed that they could increase the 
neutron star masses beyond 1.4 Mq. 

Bell and Hewish discovered radio pulsars in 1967, and they were soon identified 
as rotating neutron stars by Gold 0] . The subsequent detection of the Crab pulsar 
in the remnant of the Crab supernova, observed in China in 1054 A.D., confirmed 
the link to supernovae, and initiated the present efforts to better understand 
neutron stars. 



1.1 A Brief Overview of Observations 

Almost 1200 pulsars have been discovered by the turn of this millennium. In these 
stars the magnetic and rotational axes are misaligned, thus they emit dipole radi- 
ation in the form of radio waves that appear to pulse on and off like a lighthouse 
beacon as the pulsar beam sweeps across the Earth. The rotational energy loss 
due to dipole radiation is 

where the moment of inertia for a typical neutron star is / ~ 10^^ g cm^. Pulsars 
have magnetic fields i? of ~ 10^^ G, deduced from the observed Cl, and inde- 
pendently confirmed by cyclotron absorption lines found in X-ray spectra. Their 
periods, P = 2Tr/Q, ranging from 1.5 ms to 8.5 s, are increasing with derivatives 
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P ~ 10"^2 _ 10-21. The pulsar age is approximately given by P/2p Q; most 
pulsars are old and slowly rotating with relatively small period derivatives, except 
for a few young pulsars, e.g., those found in the Crab and Vela nebulae. 

In 1969 the Crab and the Vela pulsars were observed to "glitch", i.e. to sud- 
denly speedup with period changes AP/P of the order of 10~^ and 10~^ re- 
spectively [^. In post-glitch relaxation most of the period increase AP decays. 
These pulsars have glitched several times since then. The glitches suggest that 
the neutron stars have a solid crust containing superfluid neutrons. The inter- 
esting structure of their crust has been recently reviewed and we discuss it 
rather briefly in this report. 

The first binary of two pulsars was found by Hulse and Taylor in 1973 and they 
could determine many of its parameters including both masses, orbital period 
and period derivative, orbital distance and inclination. General relativity could 
be tested to an unprecedented accuracy by measuring the inward spiralling of 
the neutron stars in the Hulse- Taylor binary PSR 1913+16 [|ll| ]. The periastron 
advance in PSR 1913+16 is 4.2° per year as compared to 43" per century for 
Mercury, which originally was used by Einstein to test his theory of general 
relativity. Six double neutron star binaries are known so far, and neutron stars 
in all of them have masses in the range 1.36 ± O.OSMq [12]. They confirm that 
nuclear forces have a large effect on the structure of neutron stars and increase 
their maximum mass beyond 1.4 Mq. Neutron stars are estimated to have a 
binding energy of ^ 10% of their mass. Thus ~ 1.5 Mq of nuclei are needed to 
obtain a 1.35 Mq star. 

A distinct subclass of radio pulsars are millisecond pulsars with periods ^100 
ms. The fastest pulsar known has a period of 1.56 ms |1S]. The period derivatives 
of millisecond pulsars are very small corresponding to low magnetic fields ~ 
10^ — lO^^G. They are believed to be recycled pulsars, i.e. old pulsars that have 
been spun up by mass accretion whereby the magnetic fields have decayed. About 
80% of the millisecond pulsars are in binaries whereas less than 1% of normal 
radio pulsars are in binaries. About 20 - almost half of the millisecond pulsars - 
are found in binaries where the companion is either a white dwarf or a neutron 
star. 

With X-ray detectors on board satellites since the early 1970's about two hun- 
dred X-ray pulsars and bursters have been found of which the rotational period 
has been determined for about sixty. The X-ray pulsars and bursters fl^ are 
believed to be neutron stars accreting matter from high (M^IOMq) and low 
mass {M^1.2Mq) companions respectively. The X-ray pulses are attributed to 
strong accretion on the magnetic poles emitting X-rays (as northern lights). The 
observed radiation is pulsed with the rotational frequency of the accreting star. 
X-ray bursts are thermonuclear explosions of accreted matter on the surface of 
neutron stars. After accumulating hydrogen on the surface for hours, pressure 
and temperature become sufficient to trigger a runaway thermonuclear explosion 
seen as an X-ray burst that lasts a few seconds [|l^ . Masses of these stars are 
less accurately measured than for binary pulsars. We mention recent mass de- 
terminations for the X-ray pulsar Vela X-1: M = 1.87l:o;i7 Mq |15|, and the 



burster Cygnus X-2: M = 1.8 ± 0.4)Mq ||1^. They are larger than the typi- 
cal 1.36 ± O.O8M0 masses found in pulsars binaries, presumably due to accreted 
matter. 

A subclass of half a dozen anomalous X-ray pulsars has been discovered. They 
are slowly rotating, P ~ 10 sec, but rapidly slowing down. This requires huge 
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magnetic fields of S ~ 10^"^ G and they have appropriately been named "mag- 
netars" |17|. Four gamma ray repeaters discovered so far are also believed to 
be slowly rotating neutron stars. The magnetars and likely also the gamma ray 
repeaters reside inside supernova remnants. 

Recently, quasi-periodic oscillations (QPO) have been found in 12 binaries 
of neutron stars with low mass companions. If the QPO originate from the 
innermost stable orbit [18, ll^ of the accreting matter, their observed values 



imply that the accreting neutron stars have masses up to ~ 2.3Mq. In this case 
the QPO's also constrain the radii of the accreting star. 

Non-rotating and non-accreting neutron stars are virtually undetectable but 
the Hubble space telescope has observed one thermally radiating neutron star 
[pc| ]. Its surface temperature is T ~ 6 x 10^ K~ 50 eV and its distance is 
less than 120 pc from Earth. Circumstantial evidence indicate a distance of 
~ 80 pc which leads to a radius of 12-13 km for this star. In recent years much 
effort has been devoted to measuring pulsar temperatures, especially with the 
Einstein Observatory and ROSAT. Surface temperatures of a few pulsars have 
been measured, and upper limits have been set for many |21|. 

^From the human point of view supernova explosions are rare in our and neigh- 
boring galaxies. The predicted rate is 1-3 per century in our galaxy and the most 
recent one was 1987A in LMC. No neutron star associated with this explosion 
has been detected; however, 19 neutrinos were detected on earth from 1987A 
[p^ j, indicating the formation of a "proto-neutron star". It has been suggested 
by Bethe and Brown [^] that an upper limit to the mass of neutron stars can be 
obtained assuming that the remnant of SN 1987A collapsed into a black hole. 

Astrophysicists expect a large abundance of ~ 10^ neutron stars in our galaxy. 
At least as many supernova explosions, responsible for all heavier elements present 
in our Universe today, have occurred. The scarcity of neutron stars in the solar 
neighborhood may be due to production of black holes or other remnants in 
supernovae, or due to a high initial velocity (asymmetric "kick") received during 
their birth in supernovae. Recently, many neutron stars have been found far 
away from their supernova remnants; and of the ~ 1200 discovered radio pulsars 
only about ~ 10 can be associated with the 220 known supernova remnants. 
Neutron stars thrown out of the galactic plane may be detected by gravitational 
microlensing experiments [24| designed to search for dark massive objects in the 
galactic halo. 

The recent discovery of afterglow in Gamma Ray Bursters (GRB) allows de- 
termination of their very high redshifts {z > 1). They imply that GRB occur at 
enormous distances. Evidence for beaming has been observed |25], and the esti- 
mated energy output is ~ 10^^ ergs. Such enormous energies can be produced in 
neutron star mergers eventually forming black holes. From abundance of binary 
pulsars one can estimate the rate of neutron star mergers; it is compatible with 
the rate of GRB of approximately one per day. Another possible mechanism, is 
a special class of type Ic supernova (hypernovae) where cores collapse to black 
holes 

The future of neutron star observations looks bright as new windows are about 
to open. A new fleet of X- and Gamma-ray satellites have and will be launched. 
With upgraded ground based observatories and detectors for neutrinos and grav- 
itational waves 1 27] our knowledge of neutron star properties will be greatly im- 
proved. 
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1.2 Theory of Neutron Star Matter 

Neutron stars are made up of relatively cold, charge neutral matter with densities 
up to ~ 7 times the equilibrium density po = 0.16 nucleons/fm^ of charged nuclear 
matter in nuclei. The matter density is > po over most of the star, apart from the 
relatively thin crust [|T^ . The Fermi energy of neutron star matter is in excess of 
tens of MeV, and hence, at typical temperatures of ^ KeV, thermal effects are a 
minor perturbation on the gross structure of the star. 

Matter at such densities has not yet been produced in the laboratory, its prop- 
erties must be theoretically deduced from the available terrestrial data with guid- 
ance from observed neutron star properties. The quantities of interest are the 
phase and composition of cold catalyzed neutral dense matter, its energy density 
e(p) and pressure P{p), where p denotes the baryon number density. The baryon 
number is conserved in all known interactions, therefore it is convenient to find 
the composition by minimizing the total energy E^ip) per baryon, including rest 
mass contributions. This gives: 

e{p)=pEr{pl Pip) = p'^- (2) 

The equation of state (EOS) -P(e) is found by eliminating p from the above two. 

The gravitational equilibrium of a nonrotating star is described by the Tolman- 
Oppenheimer-Volkoff (TOV) § Eq: 

dP{r) _ G{e{r)+P{r)/c^){m{r)+ATTr^P{r)/c^) 
dr r2(l -2Gm(r)/rc2) ' ^' 

where G is the gravitational constant, P{r) and e(r) are the pressure and mass 
density at radius r in the star, and 

m(r) = I 47rr'^e(r')dr' , (4) 



is the mass inside r. If we neglect the general relativistic corrections of order 
the TOV Eq. reduces to the Newtonian hydrodynamic equation. The TOV Eq. 
can be easily integrated starting from the central density ec ar r = to find the 
density profile e(r). At the radius R of the star P{R) = 0, and m{R) = M is the 
mass of the star as seen from outside. The stability of the star can be deduced 
from the M(ec) as discussed in Q, and the equations for rotating stars are given 



by [28|. The effect of rotation on the structure of most observed neutron stars 



seems to be rather small, however, it could be significant at periods less than a 



millisecond |2S] 



At densities < 2 x 10~^yOo matter is believed to have the form of a lattice of 
nuclei in a relativistic degenerate electron gas |jl^ , qualitatively similar to that of 
metals. The main focus of the theory reviewed here has been on determining the 
properties and EOS of matter in the density range 2 x 10~^po < P < lOpo from 
terrestrial data. In the lower part of this range we expect to find nucleon matter 
(NM) composed of nucleons and electrons. In contrast to matter in nuclei, it has 
mostly neutrons with a small fraction of protons and equal number of electrons 
to maintain charge neutrality. The large Fermi energy, p^ ~ 100 MeV, of the 
electron gas limits the fraction of protons in NM. 

At higher densities there are several possibilities including condensation of 
negatively charged pions and kaons, occurrence of hyperons, and the transition 
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from hadronic to quark matter. All these possibilities exploit the large electron 
Fermi energy of NM, therefore only one of these, if any, may occur and lower the 
/ig- In addition, neutron star matter can have interesting mixed phase regions in 
which the mixing phases are charged but the matter is overall neutral [p!o| . 

We begin with a review of NM, and later consider the more exotic possibilities. 
In the last sections the range of neutron star structures predicted by theory is 
presented along with a comparison with the observational data. 



2 Energy-Density Functionals of Nucleon Matter 

The simplest description of nuclei is obtained within the mean field approxima- 
tion. It assumes, following the nuclear shell model, that nucleons occupy single 
particle orbitals in an average potential well produced by nuclear forces. The 
energy of the nucleus is assumed to be a functional of the orbitals occupied by 
the nucleons, and the orbitals are determined variationally as in the Hartree- 
Fock approximation. In reality the mean field approximation is not exactly valid 



for nuclei. The observed proton knockout reaction rates [30| indicate that the 
shell model orbitals are occupied with a probability of ~ 70% in the simplest 
closed shell nuclei like ^^^Pb. The differences between the real and the mean field 
wave- functions, due to correlations induced by nuclear forces, are subsumed in 



the energy functional as suggested by Kohn |31] in the context of atomic and 
molecular physics. 

The energy density of hypothetical, uniform NM at zero temperature is the 
main term in the energy functionals. The nucleon orbitals in uniform matter 
are simple plane waves, and the ground state in mean field approximation is 
obtained by filling the proton and neutron states up to their Fermi momenta 
kp^j^ = (SvrVAf)^''^, where = n,p for neutrons and protons. The energy density, 
denoted by £{pn,Pp), includes kinetic and strong interaction contributions, but 
excludes rest masses and the Coulomb interaction, which destabilizes uniform 
charged matter. The total density is denoted hy p = pn + Pp, the asymmetry of 
the matter is defined as /3 = {pn — Pp)/p, and the energy per nucleon, E{p,P), is 
given by £{p,(j)/p. 

Analysis of nuclear properties with the liquid drop model Q reveals that, 
in the absence of electromagnetic forces, the ground state of NM is symmet- 
ric (i.e. 13 = 0), has total equilibrium density pQ = 0.16 ±0.01 fm^^, and binding 
energy Eq = — 16±0.5 MeV per nucleon. The symmetry energy Egymipo) = 34ib6 
MeV, is defined as ^d'^E/dp'^ at equilibrium. The NM energy E{p,j3) can 
be expanded about its minimum value at /3 = in powers of assuming 



charge symmetry of nuclear forces. In variational |32] as well as Brueckner 
[p^ j theories the coefficients of terms with /^"-^ are estimated to be small, and 
E{p,(3) (1 — (3'^)E{p,0) + 0^E{p, 1). In this approximation the symmetry en- 
ergy is the difference between the energy of pure neutron matter and symmetric 



nuclear matter. The incompressibility Kq = 240 it 30 MeV |34] of symmetric 
nuclear matter is defined as Kq = kpd'^E/dkp at equilibrium. The energies of 
the collective breathing mode vibrations of nuclei are sensitive to Kq; however, 
in all stable nuclei the surface effects are significant. It is difficult to extract the 
density and /3-dependence of the incompressibility, and the density-dependence 
of the symmetry energy from available nuclear data. 

Analysis of elastic scattering of nucleons by nuclei shows that the nuclear mean 
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field has a dependence on the energy of the moving nucleon Q. Over a wide 
energy range this dependence is approximately linear, suggesting that nucleons 
in equilibrium nuclear matter have an effective mass m* ~ 0.7m, where m is the 
free nucleon mass. This effective mass should not be identified with the Landau 
effective mass which describes the density of single particle states in a narrow 
energy interval about the Fermi energy [ p5| , 36 1. The Landau m* in uniform 
matter is difficult to extract from nuclear data, since nucleons at the Fermi energy 
are strongly coupled to nuclear surface dynamics. Some of the phenomenological 
energy functionals are chosen to fit the observed nuclear level densities at the 
Fermi energy, while others fit the value of m*(/0o,0) obtained from the energy 
dependence of the optical model potential |^7|| . 

The nonrelativistic functionals based on Skyrme effective interactions |3^] gen- 
erally contain the following terms: 

£{Pn,Pp) = t{1 + X5P)+XiP^{1 + X2P°')+ ^ [xqTn PN + Xsp%{l + X^p") . (5) 



N=n,p 



Here tn = O.dkp ^p^ /m are the kinetic energy densities, and r = r„ + Tp. The 
parameters xi to X4 and a describe the p and P dependence of the volume integral 
of the static part of the effective interaction between nucleons in matter, while 
the X5 and xq describe effective masses produced by the momentum dependence 
of the effective interaction. In principle the values of the seven parameters in 
a typical Skyrme functional are constrained by the empirically known values of 
po, Eq, Esym{po) and K and the choice made for m*(/>o,0). However, since 
the constraints are insufficient, there are many Skyrme models of the energy 
functional. 

The simple form of the functional (Eqj5|) chosen by most Skyrme models is 
convenient, but the real functional can be much more complex. The analytic 
form of the energy density predicted by realistic models of nuclear forces, as 
discussed in the next chapter, has been studied by Ravenhall [^] . A much more 
elaborate function of the type: 



£{pn,Pp) = -P^ 
--Apn- Ppf 



Pre 



-P6P 







+ P11^J 






V P 




e-P'P) + 






g(P9P)^ 




V P 





+ J2 ^n[1 + ip3P + P5pN)e-P*P] , 

N=n,p 



(6) 



is required to reproduce the predicted £{pn, Pp) up to /3 = 1 and p ~ 1 fm~'^. This 
functional also explains the nuclear binding energies and the empirically known 
values for symmetric nuclear matter [^], however, it is unlikely that the values 
of all of its thirteen parameters can be obtained by fitting nuclear data. 

The energy of NM can be easily calculated from a covariant effective Lagrangian 



in the mean field approximation, as shown by Walecka |4C], and in the past 
decade many properties of medium and heavy nuclei have been studied with this 



approach |41|. The effective Lagrangian used in the recent work has the form: 



Jr = ^P[y^' {idf, - g^ujf, - gpf ■ p^) 



m 
1 



1 3,1 4 
-^92(7 + -^930- 



pv 



(7) 
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Hear ^, oj^ and are respectively the nucleon and uj and p vector- meson fields. 
Overhead arrows are used to denote isospin vectors, and i}^'^ = d^uj'^ — d'^iv^ etc. 
The effective scalar field a is responsible for nuclear binding, and the and cr^ 
terms are necessary to obtain the empirical incompressibility of nuclear matter 
[^ ]. The isovector p field is required to obtain the empirical symmetry energy. 
The observed values of the masses m, and rup are used, and the coupling 
constants g^^j, Qp, Qa, g2 and ^3, as well as the mass of the effective scalar 
field are adjusted to fit the nuclear data. The above Lagrangian, without the 
and a'^ terms but including pion fields and their coupling to the nucleon, is also 
used to model the two-nucleon interaction discussed in the next section. 

The relativistic mean field theory of nuclei is very elegant and often used to 
study properties of neutron star matter ji^ j . It has provided important insights 
into relativistic effects in nuclei and NM. However, the effective mean-field La- 
grangian (Eq.0) is unlikely to have a simple physical meaning. The inverse masses 
of the vector and scalar fields correspond to lengths of ~ 0.25 and 0.4 fm, which 
are much smaller than the unit radius rg = (47rp/3)^^/'^ ~ 1.2 fm for equilibrium 
nuclear matter. The naive condition for the validity of the mean field approxima- 
tion, that ro be much less than the inverse masses of the fields is totally violated 
in nuclei as well as in neutron stars. Pions are omitted from the effective La- 
grangian because they do not contribute to the energy of matter in the mean field 
approximation. Their higher order contributions are subsumed in the effective 
scalar field. Therefore, the effective mean-field Lagrangian must be interpreted 
as a relativistic generalization of Kohn's energy functional. Eq. assumes the 
simplest form necessary to fit nuclear data. A more general form, necessary to 
explain the properties of NM over the wide density-asymmetry range in neutron 
stars, can have additional fields, isovector scalar for example, density dependent 
coupling constants to take into account the changes in correlations with density, 
and field energies containing high powers of the fields, etc. 

The energy of low density neutron matter is well determined by realistic models 
of two-nucleon interaction obtained by fitting the nucleon-nucleon (NN) scatter- 
ing data. Different models and methods of calculation give very similar results 
up to p ~ /9o! beyond which three- nucleon interactions and relativistic effects, as 
well as computational difficulties may become appreciable. These energies thus 
provide a test of the ability of the Skyrme and relativistic mean field theories to 
find neutron matter properties by extrapolating data on nuclear binding energies, 
sizes, vibrations, etc. As shown in Fig.l, the neutron matter energies predicted 
by the various functionals are widely different, and not in agreement with the 
results of many-body calculations at p < pQ. It thus appears likely that the 
simple forms of effective interactions or Lagrangians used in the present mean 
field theories are inadequate to predict the properties of neutron star matter by 
extrapolating the observed nuclear properties. Nevertheless, effective mean-field 
Lagrangians have been widely used in neutron star studies due to their simplicity 
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3 Many-Body Theory of Nucleon Matter 

Many properties of nuclei and nuclear matter can be understood from the Hamil- 
tonian: 

^=-E^V2 + ^^;,,.+ Y: V,,, + .--, (8) 

i i<j i<j<k 

which includes the kinetic energy, two-nucleon interactions denoted by lower case 
V and three- nucleon interactions by capital V. The ellipsis denote neglected four 
and higher body interactions. In this section we review the present status of this 
approach and its limitations. 



3.1 Models of Two Nucleon Interaction 

Our understanding of QCD has not yet progressed enough to predict the two- 
nucleon interaction Vij ah initio. The long range part of Vij is known to be 
mediated by pions, the lightest of all the mesons, and it is denoted by the one 
pion exchange potential (OPEP) given by: 

~ 47r 3 ' 
Xij = Y^^irij) CTi ■ o-j + TT,{rij) Sij . (10) 



The pion-nucleon coupling constant f^j^^/A-K = 0.075 ±0.002 [52], and the radial 



functions associated with the spin-spin and the tensor parts are: 

YAr) = — iY{r), T^{r)=(l + - + ^\Y^{T)iT{r) , (11) 



X 



where x = nij^r. The tensor operator, Sij = 3 cTj • rijCTj ■ f jj — cTj • aj , and since 
Txir) ^ ^7r(^) in the important x^l region, the OPEP is dominated by its tensor 
part. The complete Vij is expressed as vfj + vj^j, where contains all the other, 
heavy meson, multiple meson and quark exchange parts. It has to be obtained 
along with the short range cutoffs ^y(r) and ■^r(^) in the OPEP by fitting NN 
scattering data. In boson-exchange models the is approximated by a sum of 
attractive scalar meson, and repulsive vector meson exchange potentials, while 
other models use attractive two-pion exchange and repulsive core potentials. The 
OPEP contains a (5— function term omitted from Eg . (p^) . Due to the finite size 
of nucleons this term acquires a finite range, and is difficult to separate from v^. 

In the early 1990's the Nijmegen group [^] carefully examined all the data 
on elastic NN scattering, at energies below the pion production threshold of 
~ 350 MeV, published between 1955 and 1992. They extracted 1787 proton- 
proton and 2514 proton-neutron "reliable" data, and showed that these could 
determine all NN scattering phase shifts and mixing parameters quite accurately, 
thus claiming that the experimental information on elastic NN scattering is now 
complete. Additional measurements are being carried out at several laboratories 
including the Indiana University Cyclotron Facility, and the CELCIUS facility in 
Uppsala, Sweden, to test the accuracy of the claim and improve on the quality 
of the data base. Nevertheless the Nijmegen analysis has been a major step. 

NN interaction models which fit the Nijmegen data base with a X^/Ndata ~ 1 
are called "modern". These include the Nijmegen models p9] called Nijmegen I, 
II and Reid 93, the Argonne vis |5^] (A18) and CD-Bonn In order to fit 
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both the proton-proton and neutron-proton scattering data simultaneously and 
accurately, these models include a detailed description of the electromagnetic 
interactions and terms that violate the isospin symmetry of the strong interaction 
via the differences in the masses of the charged and neutral pions, etc. 

These five models use different parameterizations of the vfj, and the Nijmegen-I 
and CD-Bonn also include nonlocalities suggested by boson-exchange representa- 
tions. Thus, like the older models, they make different predictions for the many- 
body systems. However, the differences in their predictions are much smaller 
than those between older models, and can be partly understood. 

The interaction in the spin-isospin T,S = 0,1 has the largest model depen- 
dence, which gets carried over to the estimates of the energy of symmetric nuclear 
matter (sec. 3.4). Fortunately the deuteron structure provides significant infor- 
mation on this interaction. Fig. 2 shows the deuteron wave functions obtained 
with the five modern potentials [55|. The three potentials, Reid 93, Nijmegen 
II and A18, which are local in each NN partial wave, give essentially the same 
deuteron wave function. We expect that they will give rather similar matter 
properties. The Nijmegen I and CD-Bonn potentials have momentum depen- 
dent terms associated with heavy meson exchange. These two potentials give 
larger ^Si wave functions at r < 0.8 fm, because they have softer repulsive cores 
than the local models; however, this effect is not very large as can be seen from 
Fig. 2. At r > 0.8 fm only the CD-Bonn predictions differ from the rest. This is 
because CD-Bonn has a strongly nonlocal OPEP as suggested by pseudoscalar 
pion-nucleon coupling, which suppresses the ^Di wave function. 

These results indicate that the main difference between the preset models of 
Vij is from the assumed nonlocality in OPEP. However, relativistic field theo- 
ries permit use of OPEP with different nonlocalities |^] related with the Dyson 
transformation. The three- nucleon interaction, Vijk, depends upon the choice of 
OPEP [^], and the final results obtained after including it should be indepen- 
dent of the choice. Therefore, if relativistic field theories can be used to describe 
pion exchange forces, one can use either the local OPEP in A18, or the nonlocal 
one in CD-Bonn. In this case it may be better to use the local representation 
because it is simpler, and more accurate many-body calculations can be carried 
out with it. 

Realistic models of nuclear forces, the modern as well as the older, predict 
the existence of a dense toroidal inner core in the deuteron [^]. The density 
distribution, /o(r'), of the deuteron in spin projection M = state, in the center 
of mass frame is shown in Fig. 3. Here r' = ibr/2 are the nucleon coordinates 
in the center of mass frame, and r is the internucleon distance. This density 
distribution is symmetric under rotations about the Z'-axis, and the top part 
shows its cross section in the X' — Z' plane as predicted by the A18 model. 
The bottom part shows the toroidal shape of the equi-density surface for half 
maximum density. The density peaks on a ring of diameter of ~ 1.0 fm inside 
this torus having thickness of ~ 0.8 fm. The shape is produced by constructive 
(destructive) interference between the S- and the D-wave functions shown in Fig. 2 
along the X' {Z') axis of the deuteron in the M = state. 

The dominant, static part of the Vij in the deuteron, obtained by omitting the 
terms dependent on the angular momentum, is anisotropic due to the tensor part 
of the OPEP. It strongly depends upon the angle between the unit vector r^j and 
the spin directions. The expectation value of the interaction in the 5 = 1, Ms = 
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state, "^1 Ti + it)) is shown in Fig. 4 as a function of rij for 6* = and 9 = Tr/2. 
The interaction is attractive for 9 = Tr/2 and repulsive for 9 = 0, like that between 
two magnetic dipoles. The OPEP is also shown in Fig. 4 by dashed lines. The 
NN interaction in all states except those with T,S = 1,0 is dominated by the 
OPEP at rij > 1 fm. The OPEP is weakest in the r,S' = 1,0 states. At small 
rij the repulsive core of dominates in all states. 

Most of the M = deuteron wave function has also Ms = 0. Thus it is possible 
to understand the density distribution of the deuteron from the potential shown 
in Fig. 4. The two peaks in the density shown in Fig. 3 correspond to the two- 
nucleons in the deuteron being ~ 1 fm apart at ~ 7r/2 where the potential has 
its minimum value of ^ —200 MeV in the A18 model. Other models also have 
a deep minimum at this position. The smallness of the density along the Z-axis 
is due to the repulsive potential at 9 = 0. Even though these features existed 
in the potential models of the sixties their experimental confirmation came via 
a series of measurements of the electromagnetic form factors of the deuteron up 
to momentum transfers of ~ 8 fm~^ conducted since the mid eighties at SLAG 
and Bates. These and the more recent, high precision measurements |5^, ^] 
carried out at the Jefferson Lab are in good agreement with the predictions of 
the A18 model, and verify the predicted deuteron structure beyond r ~ 0.7, or 
equivalently r' > 0.35 fm. 

The T,S = 1,0 two-proton distribution functions, believed to be similar to 
the two-neutron distribution functions due to isospin symmetry, are predicted 
1 58 1 to have a dip at r ~ and a peak at r ~ 1 fm, due to the repulsive core 
and the minima of v^n respectively. The experimental information on these 
spherically symmetric distribution functions is less direct. It comes from sums 
of longitudinal response functions of light nuclei [^]. The sums observed in '^He 
and ^He are in fair agreement with theory, and show evidence of the predicted 
structure, however, the relativistic and other corrections to the observed sums 
are significant. 

The observed deuteron form factors and, to a lesser extent, the sums of longi- 
tudinal response of light nuclei indicate that the modern two-nucleon potentials 
and the wave functions they predict have validity at internucleon distances larger 
than ~ 0.7 fm. This may appear surprising because the rms charge radius of the 
proton is known to be 0.8 fm. However, the nucleons seem to have a small and 
dense core. The charge form factor of the proton, as well as the magnetic form fac- 
tors of the proton and neutron are well approximated by the dipole (l + (/^/9o)~^ 
with qo = 840 MeV/c. Inverting this form factor gives the proton charge density 
as p^fi{r) = 3.3 Q-^/'^-'^^f"^ fm~^. The charge densities of two protons, one fm 
apart are shown in Fig. 5. It should be noted that the charge densities shown 
can have corrections due to the neglect of relativistic effects in inverting the form 
factor at r^l/m ~ 0.21 fm, and that recent measurements of the proton form 
factor [S3 1 show deviations from the dipole form at momenta ^ 4 fm^^. They 
suggest that proton charge density flattens out at r < 0.3 fm. Nevertheless Fig. 
5 indicates that nucleons a fm apart can retain their identities. The NN inter- 
action includes the change in the energy due to their overlap, and it has minima 
near r ~ 1 fm. In absence of the quantum kinetic energy term (— V^/2m) in the 
Hamiltonian (Eq.^) the deuteron will shrink to a ring of radius ~ 0.5 fm, and the 
equilibrium density of nuclear matter will be ~ 1 fm~^. The density of matter 
in most neutron stars is less than that. 
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3.2 Models of Three Nucleon Interaction 

All realistic models of Vij, the modern and the older, underbind the triton and 
other light nuclei and predict too high equilibrium density for symmetric nuclear 
matter. In both cases the deviation from experiment is not too large, particularly 
when compared with the expectation values of Vij. For example, the expectation 
value {vi<^j) in is about ~ —50 MeV, while the underbinding is by < 1 MeV. 
It is likely that these differences are due to three-nucleon interactions expected 
and predicted since the fifties 

High precision modeling of the two-nucleon interaction is possible because scat- 
tering cross sections can be easily and exactly calculated from the assumed u,-,-, 
and a complete set of ~ 4000 cross sections has been measured. Such an ap- 
proach is not practical at present for the three-nucleon interaction Vijk- In prin- 
ciple deuteron-nucleon scattering and reactions can be used to study Vijk in the 
isospin T = 1/2 state. However, this scattering is dominated by the two-nucleon 



interaction [65|, and very high precision data is necessary to extract the effects 
of Vijk, as is being attempted in new experiments at the Indiana University cy- 
clotron, focused on spin observables. 

The large ( Vi^j) is cancelled to a large extent by the kinetic energy in nuclear 
binding energies. Thus (Vi<j<fc) in nuclei is expected to be of the order of 10 
% of their binding energy. Hence we can construct realistic models of Vijk by 
fitting binding energies of light nuclei, which can now be calculated with an ac- 
curacy of the order of 1 % using Greens Function Monte Carlo (GFMC) methods 
1^, ^ and the estimated equilibrium properties of nuclear matter. As will be 
discussed in sec. 3.7, it is not yet possible to calculate the EOS of symmetric 
nuclear matter with comparable accuracy. The models depend upon the Vij used 
in the Hamiltonian. This is inevitable, because unitary transformations make 
correlated changes in Vij and Vijk Only the combinations of Vij and Vijk in 
the Hamiltonian (|8|) are meaningful. 

The information contained in nuclear binding energies and equilibrium prop- 
erties of nuclear matter is limited. Therefore the realistic models of Vijk I'ely 
on theory to a much larger extent than the models of Vij, and contain very few 
parameters. The Urbana models of Vijk contain two isoscalar terms: 

V^,k = V;^l + V^, . (12) 

The first term represents the Fujita-Miyazawa |64] two-pion exchange interaction: 

A27r ■ ' - 1 , 



eye 



^ijk = 2^ ^27r [ {Tj • Tj, Tj • Tfe} {Xij , Xik} + -[Vi ■ Tj,Ti ■ Tk][Xij,Xik 



(13) 

with strength denoted by A2n- The functions Tj^irij) and y^(rjj) in Xij, Eq. (|lO|) , 
are taken from the A18 model of Vij. This interaction is due to the pion exchanged 
by nucleons j and k being scattered by the nucleon i via the A resonance in vr- 
N scattering. In classical terms it is due to the polarization of the quark spins 
in nucleon i, due to the pion field of j (k), interacting with k (j); and it is 
similar to the three-body earth-moon-satellite gravitational interaction due to 
the polarization of the ocean water on earth by the moon's gravity. 
The Vj^f^ is purely phenomenological, and has the form: 

Vifk = UoJ2T'{rij)T^{nk). (14) 

eye 
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Table 1: Results of Quantum Monte Carlo Calculations in MeV 
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It was meant to represent the modification of the two-pion exchange part of 
Vij by other particles in matter, however, ~ 40 % of it is due to relativistic 
effects discussed in the next subsection. The parameters of the present model 
U-IX, ^27r = -0.0293 MeV and Uq = 0.0048 MeV, have been determined from 
exact GFMC calculations of '^H and approximate variational calculations of the 
equilibrium density of nuclear matter with the A18 NN interaction 



The results of essentially exact GFMC calculations |6^, 67] with the A18 and 
U-IX interactions are shown in Table-|I[ We note that the better known pion 
exchange parts of these interactions give the largest contributions, but the con- 
tributions of the phenomenological parts, and V^jj^, are significant. The column 
AEexpt gives the difference between the experimental energies and the calculated, 
while AEvMC is that between the variational Monte Carlo (VMC) upper bounds 
and the exact GFMC energies. 

Table-|l| shows that the U-IX interaction underbinds A = 8 nuclei, and since 
'^He is more underbound than ^Be, it misrepresents the isospin dependence of 



Vijk- The new Illinois models of Vijk [69| resolve this problem by including the 



leading three-pion exchange term, V^'^ , that is attractive in triplets having isospin 
T = 3/2, but has little effect on the T = 1/2 triplets in and ^He. A much 
improved fit, with errors < 2%, to the observed energies is obtained as shown in 
Fig. 6. The three parameters in Pieper's model IL-2R are the strengths of the V^^, 
V^^j, and V^J.. Calculations of NM properties with this more accurate Vijk are in 
progress; the preliminary results are similar to those with U-IX since the V^'^ is 
much weaker than the V'^'^ and V^. For example, the expectation values of V^'^., 

V^,^ and V^l of IL-2R model, in ^Be (^He) are respectively 38 (-27), +19 

(+14) and —2 (—5) MeV. In the following sections we will review the properties 
of neutron star matter calculated with A18 and U-IX interactions. 

In principle there can be four-nucleon interactions (FNI) neglected in the 
Hamiltonian (Eq.|8|). It seems that they are very weak in nuclei. All the models 
of Vijk studied so far reproduce the energy of ^He with an error < 0.5%, after 
fitting the observed energy of ^H. Since this error is close to the accuracy of the 
^He calculation, there is no indication of FNI in that nucleus. The IL-2R model 
also gives the experimental energies oi A = 8 nuclei within ~ 1 MeV. In ^Be, for 
example, the expectation values of Vij and Vijk are respectively —308 and —21 
MeV respectively. By comparison with experiment we estimate that the possi- 
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ble contribution of FNI in this nucleus is < 1 MeV. The Vijk presumably has 
additional smaller terms neglected in IL-2R, but it is difficult to determine their 
strengths from the nuclear spectra that can be calculated accurately from bare 
forces at present. 



3.3 Relativistic Boost Interaction 

In all models, the NN scattering data is reduced to the center of mass frame and 
fitted using phase shifts calculated from the NN interaction, Vij, in that frame. 
The Vij obtained by this procedure describes the interaction between nucleons 
having total momentum Pjj = Pi + Pj = 0. In general, the interaction between 
particles depends on their momenta pj and pj. For example the electromagnetic 
interaction |^0| between two particles of mass m and charge Q, contains a term 



with the factor — pj •pj/2m^ = pfj/2m'^ — P^j/8m?, where Pij = {pi — pj)/2 is the 
relative momentum. The terms containing pij are included in the momentum- 
dependent parts of Vij, while the boost interaction 6v{Pij) contains parts depen- 
dent on the total P jj . Even though contributions of the boost interaction to the 
binding energy of SNM and were estimated by Coester and coworkers in the 



seventies and eighties |7l|, ^] , they were first included in studies of dense matter 
rather recently [73|. Walecka's relativistic mean field theory naturally contains 
the boost interactions ff^ . 

Following the work of Krajcik and Foldy, Friar [^] obtained the following 
equation relating the boost interaction of order to the interaction in the center 
of mass frame: 

The general validity of this equation in relativistic mechanics and field theory 



was recently demonstrated [74|. Including boost interaction, the nonrelativistic 



Hamiltonian assumes the form: 

^* = E ^ + E(%- + '^^(p^.)) + E + ■ • • ' (16) 

where the ellipsis denotes the three-body boost, and four and higher body inter- 
actions. This H* contains all terms quadratic in the particle velocities, and is 
therefore suitable for complete studies in the nonrelativistic limit. 

Studies of light nuclei using the VMC method |[76| , 77 1 find that the contribution 
of the two-body boost interaction to the energy is repulsive, with a magnitude 
which is ~ 37% of that of Vjji^ in the UIX model. The boost interaction thus ac- 
counts for a significant part of the V-^j^ in Hamiltonians which fit nuclear energies 
neglecting 6v. The V^*^ in Eq.(^6|) has a of appropriately smaller strength 
than that in the Vij^ in Hamiltonian H given by Eq.(^. 

We should expect additional relativistic corrections to the Hamiltonian (16). 



However, when nonrelativistic potentials are fit to the experimental data, rela- 
tivistic effects present in the data are automatically buried in these potentials. In 
order to study the magnitude of a chosen relativistic correction, such as that due 
to the approximation of the kinetic energy or the boost interaction, or the non- 
localities of OPEP, it is necessary to refit the same data set, fitted to obtain the 
nonrelativistic Hamiltonian, and then study the differences. Such comparisons 
Jr^] indicate that the relativistic corrections associated with kinetic energies and 
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nonlocalities of OPEP are small, whereas the boost corrections are significant. 
This is not surprising since the boost interaction was totally omitted from the 
conventional nonrelativistic nuclear Hamiltonian (Eqj^. 

3.4 Brueckner Calculations of Nucleon Matter 

Calculating the properties of matter from the interaction Vij between pairs of its 
constituents is a well known problem in many-body theory. It is particularly chal- 
lenging for NM due to the strong spin-isospin dependence of Vij. In the method 
developed by Brueckner, Bethe and Goldstone the perturbation expansion of the 
energy of NM is cast into a series ordered according to the number of indepen- 
dent hole lines (HL). This method has been used extensively since the sixties [^] 
to study symmetric nuclear matter (SNM) and nuclei. It was also used then to 
predict properties of pure neutron matter (PNM) ||7^ soon after the discovery 
of pulsars. The convergence of the expansion depends upon the choice of the 
single-particle energies in the assumed unperturbed Hamiltonian. For the hole 
states with momenta less than kp, they are chosen self consistently as suggested 



by Brueckner and Gammel |80| via the Brueckner-Hartree-Fock (BHF) proce- 
dure. The older calculations by Day ^l|] used kinetic energies for particle states. 
Since this leads to a discontinuity in the unperturbed single-particle energies at 



kp, Day's choice is called "discontinuous". In 1976 the Liege group |Q advo- 
cated the "continuous" choice by extending the definition of BHF single-particle 
energies to particle states. If all the higher order terms of the HL expansion are 
computed the final results should be independent of the choice. 

Detailed calculations of SNM have been carried out with the older Argonne 
vi4 interaction by the Catania group using both choices. The results of the 
lowest order 2-HL calculation depend significantly on the choice, however, those 
including 2+3-HL terms are almost the same for the two choices. They find that 
the 2-HL results with the continuous (discontinuous) choice are ~ 15 % below (30 
% above) those of the 2+3. Day's calculations |^l|, with the same interaction, 
include additional 4-HL terms, and give energies below the 2+3 Catania results, 
closer to the 2-HL continuous. For example, the energies of SNM at p = 0.28 fm~^ 
are —11.3, —16.1, —18.3 and —17.8+1.3 for the Catania 2-HL discontinuous, 2+3- 
HL, 2-HL continuous, and Day's calculations respectively, while the variational 
upper bound at this density is —16.2 + 0.4 |^l| with the methods described in 
the next subsection. The size of the error due to truncation of the expansion is 
estimated in Day's and variational calculations. 

A significant advantage of Brueckner's method is that it can be easily applied 
to local as well as nonlocal modern interactions. The present calculations, called 
lowest-order BHF (LOBHF), include only the 2-HL terms and use the continuous 



choice [B3|. Their results for PNM and SNM are shown in Fig. 7; those for matter 
with intermediate values of proton fraction can be estimated from these by inter- 
polating with as metioned in sec. 2. These results provide an estimate of the 
uncertainty in the predicted matter energy due to that in the NN interaction. 

The LOBHF energies for neutron matter are essentially model independent up 
to ~ 0.3 fm~^ (Fig. 7); at higher densities they deviate partly for the following 
reason. All the models fit the NN scattering data up to 350 MeV lab energy, i. e. 
up to maximum relative momenta k = 2.05 fm~^. The maximum value of the 
relative momentum of two hole states in matter is kp, and the density of PNM 
(SNM) at kp = 2.05 is 0.29 (0.58) fm'^. In PNM at p > 0.29 the interactions 
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are being used at relative momenta larger than in the fitted data. Baldo et al 
p^ ] have shown that the phase shifts predicted by the five modern potentials 
vary from 8 (A18) to 19 (Nijmegen II) degrees at /c = 3 fm~^. The average 
value, krms = 1/3/10 kp-, is smaller than kp, and exceeds 2.05 fm~^ at much 
larger densities of 1.77 (3.54) fm~^ of PNM (SNM). Presumably, this helps to 
keep down the model dependence. 

The LOBHF energies of SNM have a larger model dependence starting at lower 
densities (Fig. 7). Here the main cause seems to be the assumed nonlocality of 
Vij discussed in sec. 3.1. The local interactions, Nijmegen II, Ried 93 and A18 
give similar results, while the most nonlocal CD-Bonn gives the lowest energies. 
The predicted values of equilibrium pQ and Eq of SNM are respectively 0.31, 0.27, 
0.28, 0.27 and 0.37 fm'^ and -20.3, -17.6, -18.7, -18.1 and -22.9 MeV with 
Nijmegen I, II, Reid93, A18 and CD Bonn interactions; while the empirical values 
are 0.16 fm~^ and —16 MeV. Obviously, the empirical properties of SNM can not 
be obtained by approximating the nuclear interaction energy by J2'^i<j- TNI, 
added to the Hamiltonian to obtain the observed properties, naturally depend 
upon the choice of Vij. For example, those to be added to CD Bonn have to 
have stronger repulsive parts, which dominate at large p. The combinations of 
Vij+Vijk, constrained with experimental data, will have smaller model dependence 
than seen in Fig. 7. 

Baldo, Bombaci and Burgio ||8^ have carried out LOBHF calculations with 
the older Paris and Argonne interactions including Urbana TNI. Like the modern 
models, these older models give too large po (~ 0.26 fm~^) and (-^ —18 MeV) 
in LOBHF without TNI. By averaging over the position of the third nucleon the 
TNI is expressed as a density dependent NN interaction to be added to the Vij. 
The parameters and Uq were chosen to get closer to the empirical values of 
Po and £"0; their values of —0.0329 and 0.00361 MeV are not too far from those 
of U-IX (-0.0293 and 0.0048). With the Paris + Urbana model they obtain 
Po = 0.176 fm-3, Eo = -16.0 MeV, K = 281 MeV, and E,ym = 33 MeV. These 
values, as well as those obtained using the Argonne interaction instead of Paris, 
compare rather well with the empirical values given in sec. 2. 

Relativistic effects are included in the LOBHF calculations via Dirac-Brueckner 



approximation suggested by Celenza and Shakin |86]. The calculations include 
contributions of the boost interactions as well as TNI |87], and many nucleon 
interactions generated via the Z-graphs containing anti-nucleon lines. However, 
they do not contain contributions of the Fujita-Miyazawa and other TNI due to 
internal structure of nucleons. Results have been reported by Brockmann and 
Machleidt for the older Bonn meson exchange NN interaction models; those 
with the Bonn-A model come close to reproducing the empirical properties of 
SNM at Po- 

3. 5 Variational Calculations of Nucleon Matter 

Variational calculations of NM with realistic interactions have been carried out 



since 1970 [89|. The present calculations |68, 73| use variational wave- functions, 
consisting of a symmetrized product of pair correlation operators, Fij, op- 
erating on the Fermi gas wave-function. In PNM, the Fij include four terms 
generating spatial, (Xj • crj, tensor and spin orbit correlations. The SNM Fij have 
eight terms; the additional four have Tj • tj factors. 

This wave-function is clearly too simple to accurately describe the ground state 
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of nuclear matter. Monte Carlo studies of few-body nuclei use additional three- 
body correlations induced by both Vij and Vijk, in the variational wave- function; 
they reduce the energy of ^^O by ~ 1 MeV/nucleon. The results shown in Ta- 
ble.|| indicate that the VMC energies oi A = 8 nuclei, obtained after includ- 
ing three-body correlation operators, are above the exact GFMC values by ~ 1 
MeV/nucleon. From these we estimate that the present may underbind SNM 
by a few MeV. In contrast, the three-body correlations have a smaller effect on 
the energy of pure neutron drops |^^. The variational energy of a drop with 
eight neutrons, calculated with the simple is greater than the exact value 
by ~ 0.5 MeV/nucleon. Thus the variational energies are relatively more accu- 
rate for PNM than for SNM. This is as expected, since SNM has much stronger 
tensor correlations. Despite the aforementioned shortcomings, the simple 
having only pair correlation operators describes the gross features of the nuclear 
wave-function rather well. For example, the spin-isospin dependent two-nucleon 
distribution functions calculated in this approximation are close to the exact 
distribution functions [|6^ ]. 

The correlation operators Fij are determined from Euler-Lagrange equations 
that minimize the two-body cluster contribution of an interaction Vij — Xij. The 
interaction Vij is related to the Vij via a parameter a meant to simulate the 
quenching of the spin-isospin interaction between particles i and j, via their in- 
teraction with other particles in matter. The operator Xij simulates screening 
effects in matter; it is determined from the ranges dt of tensor correlations, and 
dc of all the other correlations. The thus depend on three variational pa- 
rameters: a, dc and d^, determined by minimizing the energy. Two additional 
parameters are used in Ref. |91] to further lower the variational upper bound by 
small amounts. 

The energy expectation value is evaluated using cluster expansion. The one- 
body term is just the Fermi gas kinetic energy, and the large two-body (2B) 
term, analogous to the interaction energy in LOBHF, is calculated exactly. The 
most important of the many-body (MB) cluster contributions are summed using 
Fermi-hypernetted chain (FHNC) and single operator chain equations |68, 
and constrains are imposed to satisfy the fundamental identities of pair distri- 
bution functions. The kinetic energy can be calculated using different expres- 
sions related by integration by parts. If all MB contributions are calculated, 
these expressions yield the same result. However, they yield different results 
when only selected parts of the MB clusters are summed using chain equations. 
Studies of atomic helium liquids with FHNC summation methods find the exact 
result to be between the energies obtained using the Jackson-Feenberg (JF) and 
Pandharipande-Bethe (PB) expressions. The average of these two is used as the 
result with half the difference as an estimate of the error. 

More general pair correlations can be calculated by separately minimizing the 
two-body cluster contribution to each partial wave, specified by /, S, J and the 
relative momentum k These correlations f{l,S, J, k) depend on all the quan- 
tum numbers, and yield a lower 2B energy than the Fij operator with the same 
a, dc and dt- The MB contributions cannot be easily calculated with the general 
/(/, S*, J, k), however. The differences between optimum /(/, S, J, k) and Fij can 
be included via the second order two-particle, two- hole contribution, in 
correlated basis perturbation theory In recent calculations |65, 73] the 

Ai?2 is approximated by the difference 5E2b between the 2B cluster energies cal- 
culated using f{l, S, J, k) and Fij. However, the values of a, dc, dt are determined 
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by minimizing the energy calculated from the Fij. 

The variational [^] and LOBHF energies obtained with the A18 interaction 
are compared in Fig. 8. At densities below 0.6 fm"^ there is fairly close agreement 
between them, however, we expect the true results to be a few MeV below the 
variational upper bound. At higher densities the SNM LOBHF energy is signif- 
icantly below the variational bound. The convergence of the HL expansion is 
expected to deteriorate at higher densities, and it may be the cause of the large 
difference. 



3.6 Neutral Pion Condensation 

Variational calculations with the Hamiltonian ( [l^ ) with AlS+dv+VYX* interac- 
tions, indicate the occurrence of a phase transition in both PNM and SNM; the 
energies of the two phases are shown in Fig. 9. The tensor correlations have a 
longer range in the higher density phase (HDP) than in the low density phase 
(LDP). Detailed studies of the pion fields in the two phases j68U indicate that 
the HDP has a large enhancement of virtual neutral pions with momenta ~ 1.5 
fm~^, and therefore this transition is believed to be due to neutral pion conden- 
sation in matter. Although the effect of this type of transition on the EOS is 
relatively small, it can have important consequences for the cooling and evolution 



of neutron stars [21|. 



Since the pioneering work of Migdal |M| and of Sawyer and Scallapino |95|, 
many investigators have used effective interactions to study the possibility of 
pion condensation in SNM and PNM. These efforts were recently reviewed by 
Kunihiro et al [^]. Neutral pion condensation occurs when PNM (SNM) be- 
comes unstable towards the development of a spin (spin-isospin) density wave as 
discussed in sec. 4.2. The variational wave functions used in the Urbana calcu- 
lations 1 68 1 are not adequate to describe the long range order expected with vr'^- 



condensation; better calculations may be possible with the quantum Monte Carlo 
method described in the next subsection. NM with spin or spin-isospin density 
wave naturally has a pion field of the same wavelength; it causes ^ A transi- 
tions, which help lower the energy of this phase. The effects of the A-resonance 
are absorbed into the interactions in the Hamiltonian (|T^). In particular, the 
Fujita-Miyazawa V^^^ has a large effect on condensation. Without it there is 
no condensation predicted in SNM, while that in PNM occurs at a higher density 
of ~ 0.5 fm~^. 

This transition was found by Wiringa, Fiks and Fabrocini to occur with the 
older Argonne fi4 for PNM but not for SNM |91], while it does not occur in either 



PNM or SNM with the Urbana vi4 interaction of 1981. In Migdal's approach 
the transition of SNM to the pion condensed phase is inhibited by a positive, 
short-range <Tj • crjTi ■ Tj NN interaction whose strength is represented by the 
Landau parameter g' . In the case of PNM the sum of the <Tj • crj and cr, • CTjTi ■ Tj 
interactions occurs since Ti ■ tj = 1. The contact part of OPEP gives a negative 
<Tj • (TjTi ■ Tj potential at small r, and thus favors pion condensation. The cr, • 
CTjTi ■ Tj part of the modern A18 NN interaction does become negative at small 
r, however, it is positive in the older Urbana-Argonne models which do not fit 
the NN scattering data as well as the A18. 

The energies of /3-stable LDP and HDP phases have been calculated by inter- 
polation between PNM and SNM |7|]. The LDP phase with p = 0.204 fm'^ 
and proton fraction Xp = 0.073 is found to be in equilibrium with the HDP at 
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p = 0.237 and Xp = 0.057. In between there will be mixed phase regions as 
discussed in sec. 5. 

3.7 Quantum Monte Carlo Calculations 

The variational energy of the ground state of SNM, calculated with the Al8+5v+UIX* 
model is ~ —12 MeV, against the empirical value of —16 MeV. As mentioned ear- 
lier, this difference is believed to be mostly due to the inadequacy of the present 
variational wave functions. Accurate calculations of the energies of many-body 
systems in which the interactions do not depend upon the spins of the parti- 
cles have been carried out using quantum Monte Carlo (QMC) methods [p7|]. 
In these systems one can work with a wave function <I*(R) that depends only 
upon the positions of all the particles represented by the configuration vector 
R = ri, r2, ...fyi. A VMC calculation is used to obtain a good approximation 
^>^(R), and in the GFMC method one operates on <I>^(R) with the imaginary 
time evolution operator exp{—[H — Eq\t) to project out the exact ground state 
$o(R). 

The main difficulty in applying QMC methods to nuclear problems is that 
nuclear forces change spins and isospins of the interacting nucleons, and thus 
nuclear wave functions contain superpositions of all possible spin-isospin states 
of A-nucleons. Their number, ~ 2^ A\/[Z\{A — Z)l], increases very rapidly with 
A. For this reason exact QMC calculations have been carried out only for nuclei 
having up to 8 nucleons and attempts to calculate A = 9 nuclei are in 



progress. Carlson ||9g] has also calculated the ground states of 14 neutrons in a 
periodic box. Calculations of pure neutron systems are simpler because there is 
only one isospin state, and all Ti ■ tj can be replaced by the unit operator. 

In addition, GFMC calculations of Fermion systems suffer from the "Fermion 
sign problem" . The real wave functions of simple Fermi systems have nodal sur- 
faces because the ^*(R) must equal — $(R') when the configurations R and R' are 
related by the exchange of a pair of particles. GFMC configurations which diffuse 
across nodal surfaces, as the system evolves in imaginary time, increase the vari- 
ance of the calculated quantity, making unconstrained propagation impractical 
for large systems. In the fixed node method ||9^ for simple systems this growth 
of variance is eliminated by restricting the configurations to domains enclosed by 
the nodal surfaces of the variational wave function. Such calculations generally 
have an error due to imperfections in that structure, but it is much smaller than 
that in variational. 

A similar problem comes in nuclear GFMC with the additional complexity 
due to nuclear wave functions having many spin-isospin components, each with 
a different nodal structure. The growth of the variance is tolerable for A < 
7, but when A > 8 it is necessary to use constrained path methods to 
control the variance. The constraint can be removed at large r to test if it 
influenced the calculated energy significantly. It appears that calculations with 
~ 2% accuracy in the binding energy are possible in this way for systems having 
up to 14 neutrons. 



Auxiliary field diffusion Monte Carlo (AFDMC) |10C] seems to be the long- 
sought breakthrough needed to eliminate the exponential (2^) growth of spin 
states in GFMC calculations of neutron matter. This method combines two major 
themes in QMC. Auxiliary- field methods are used in the shell model Monte Carlo 



calculations |101], and several condensed matter systems in which the continuous 
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spatial degrees of freedom have been eliminated, while diffusion Monte Carlo is 
another name for GFMC. 

In the approach developed by Schmidt and Fantoni ]100| ], the spatial parts, 
(i.e. kinetic energy and spin- independent interactions), of the Hamiltonian are 
propagated as in GFMC and the spin-dependent interactions between neutrons 
are replaced by interactions of neutrons with auxiliary fields. Integrating over the 
auxiliary fields reproduces the original spin-dependent interaction. In addition, a 
constraint analogous to the fixed-node approximation in GFMC is introduced, by 
requiring that the real part of the overlap with a trial function remains positive. 

More recently, Schmidt and Fantoni [ll02[| have carried out calculations with a 
realistic Hamiltonian consisting of Argonne NN interaction used in the GFMC 
calculations |66|, and the UIX TNI. The contains the main parts of A18, and 
the difference between the two is treated perturbatively. Results of calculations 
with 38 neutrons in a periodic box with finite size corrections have been ob- 
tained. They are ~ 5% below the variational energies obtained with the methods 
described above. For example, the AFDMC and variational energies for Argonne 
v'g and UIX interactions are 21.8 ±0.1 (65.5 ±0.1) and 23.2 (68.6) MeV per neu- 
tron at p = 0.2 (0.4) fm-3. 

The trial functions used to constrain the present AFDMC calculations are 
rather simple without any spin correlations. In contrast it is possible to use 
more accurate variational wave functions with spin correlations to constrain the 
GFMC calculations. Carlson [|8| has compared AFDMC and GFMC results for 
14 neutrons in a periodic box. At ,9 = 0.15 fm~3 the GFMC energy (220 ± 
1) is about 7 % below the AFDMC result of 236.4 ± 1.5 MeV. From these we 
conclude that the variational PNM energies given in the last section may be ~ 12 
% above the exact values for the A18 and UIX interactions. The error in SNM 
E{p) is probably twice as large. As mentioned in the conclusions (sect. 7), an 
overestimation of the E{p) of neutron matter by 12 % has a rather small effect 
on the predictions of neutron star properties. 

The AFDMC is more accurate than the present variational method, and it is 
also more versatile. For example, it can be used to study matter with long range 
spin-isospin order induced by vr'^-condensation discussed in sect 4.2. 



4 Hadronic and Quark Matter 

It is likely that at high densities more general form of matter containing hadrons 
besides the nucleons, called hadronic matter (HM), has lower energy. The pos- 
sibilities are that it contains negatively charged mesons like pions or kaons, or 
other hyperons such as T,~ or A. Finally, it is expected that at a high enough 
density there will be a transition to quark matter (QM) in which the quarks are 
not clustered into nucleons or hadrons. 

The interactions between hyperons and nucleons, and between kaons and nu- 
cleons are not as well known as those between nucleons, and the energy of quark 
matter is difficult to calculate realistically. Therefore the transition densities from 
NM to HM or QM are rather difficult to calculate reliably. We review the recent 
estimates. 
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4-1 Kaon Condensation 



Kaon condensation in dense matter was suggested by Kaplan and Nelson |103], 



and has been discussed in many recent publications |104, 105]. Due to the attrac- 
tion between K~ and nucleons the kaon energy decreases with increasing density, 
and eventually if it drops below the electron chemical potential in NM, a Bose 
condensate of K~ will appear. The key quantities of interest are the electron and 
kaon chemical potentials in NM. The former is obtained from the /3-equilibrium 
condition, fi^ = fJ-n — fJ-p relating electron, neutron and proton chemical poten- 
tials. The He obtained from the A18+6v+\]TK* interactions with variational 
calculations is shown in Fig. 10. 

In neutron matter at very low densities, the interparticle spacing is much larger 
than the range of the K~n interaction, and the kaon interacts many times with 
the same nucleon before it encounters and interacts with another nucleon. There- 
fore one can use the scattering length, a^-„, as the "effective" kaon-nucleon in- 
teraction. In this low density limit the kaon energy deviates from its rest mass 



by the Lenz potential, and is given by [IOC] 



27r , , 

(^Lenz = 'mK^ aR-n P, (17) 

where niR = rnxnin/ (itt-k + "in) is the kaon-neutron reduced mass. The scatter- 
ing length extracted from data is (—0.37 — iO.57) fm; its imaginary part is due 
to the open A — vr" channel in vacuum. In the density region of interest to kaon 
condensation the kaon energy is too small for the K~ -|- n ^ A -|- vr" reaction to 
occur. Using effective lagrangians based on chiral perturbation theory Brown et. 



al. [104 1 estimate a^-^ to be —0.41 fm in absence of reaction channels. 



There are two corrections to the Lenz energy at small densities. Including 



these, the kaon energy to obtains the form |105, IOC ]: 



/ 2mK \ / 1 \ 27r 
UJ = mK+ [ ■ — Oi^-n P , (18) 

where is the inverse correlation length. The first correction factor is a rela- 
tivistic effect obtained from the Klein-Gordon equation, while the second factor 
is from the theory developed by Ericson and Ericson for propagation of mesons 
in nuclear matter. The relativistic correction decreases the kaon energy, while 
the correlation correction increases it. 

As the density increases further, and the interparticle spacing becomes of the 
order of the range of the interaction, the kaon will simultaneously interact with 
two or more nucleons and the Lenz approximation will break down. In the high 
density limit the kaon energy deviates from its rest mass by the Hartree potential: 

i^Hartree = rUK + P j V/^-„(r)d^r , (19) 



where is the K~n interaction potential. As shown in Ref. [ [L07| , the 

Hartree potential is considerably less attractive than the Lenz potential and thus 
uj Hartree > ^Lenz- The transition from the Lenz to Hartree limits has been re- 
cently studied [106| with a variety of methods including exact calculations for 



simple cubic crystal model of neutron matter. For reasonable interaction range 
the transition begins at very low densities (< O.lpo)) and the Hartree limit is 
essentially reached by ~ 3po. There are no relativistic corrections to the Hartree 
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energy of kaons condensed in the state with zero momentum ]106| provided the 
K~n interaction is dominated by the Weinberg- Tomozawa vector potential. 

The typical resent results for kaon energy in neutron star matter are shown in 
Fig. 10. The top solid line is obtained with a Wigner-Seitz [|106(| calculation for 
pure neutron matter; it is exact in both the low and high density limits, and gives 
essentially the Hartree energy at p > 3po- The next curve shows the estimated 
Hartree results for NM containing ~ 15 % protons, below that is the tVLenz for 
neutron matter, while the lowest curves uses K~N scattering amplitude calcu- 
lated in matter including relativistic, correlation and proton fraction corrections 



[105|. The recent [73| estimate for the central density of a 2.0 Mq star is ~ bpo, 
and the results shown in Fig. 10 indicate that kaon condensation is unlikely at 
densities lower than that. The heaviest stars (M = 2.2Mq) made up of NM are 
predicted to have central densities of ~ 7po, and the possibility of kaon conden- 
sation in their cores can not be ruled out. However, presence of tt~ , or S~ or 
quark drops will decrease pe as discussed in the following subsections and may 
make kaon condensation unlikely in even the most massive stars. 

4-2 Charged Pion Condensation 

Negatively charged 7r~-mesons will condense in matter when their chemical po- 



tential becomes lower than as suggested by Migdal and by Sawyer and 



Scallapino |Q. In the seventies and eighties this possibility was studied by many 



researchers. Their work has been reviewed in Ref. |96|, and we discuss it rather 
briefly. 

Fig. 10 suggests that, in absence of interactions vr" with zero momentum will 
condense in matter at a rather low density of ~ 1.5po when /ie exceeds their rest 
mass of 139 MeV. However, the 7r~n S-wave interaction is repulsive, and raises 
the energy of zero momentum pions sufficiently above the estimated pe- The 
recently found []108(| , deeply bound 7r~-nucleus atomic states, are influenced by 
this repulsion. The K~n S-wave interaction, on the other hand is attractive, 
leading to the possibility of kaon condensation discussed in the last section. 

The energy of vr", having momenta of the order of 2 fm~^, ~ 400 MeV with- 
out interactions, is reduced in matter by the 7r~n P-wave interaction due to the 
coupling of the vr^ to p-n and A^-A particle-hole states. Their energies and cou- 
plings to 7r~ are calculated with effective forces described with Landau parameters 
whose density dependence is not well established. With plausible values for the 
Landau parameters a second order transition with vr^ condensation is predicted 
at ^ 2pq; however it is not expected to have a large effect on the EOS [96|. 



The Japan group also predicts a first order transition to an interesting phase 
with both vr" and vr*^ condensation at a density of ~ 4.5po- It has a significant 
effect on the EOS. Matter in this phase has spin aligned layers as discussed 
earlier in sec. 3.6 and illustrated in Fig. 11. The condensed vr^-mesons have 
momenta in the Z-direction, perpendicular to the layers, while that of vr" is in 
the X-direction in the plane of the layers. This way, when a proton absorbs a vr^ 
form the condensate and becomes a neutron its spin direction also gets flipped 
maintaining the attractive interaction with the vr" field. It is a challenge to 
calculate the energy of matter in this interesting phase from bare nuclear forces. 

The Illinois group, working with bare nuclear forces, predicts the first order 
vr'^-condensation at a much lower density of ~ l.Spo (sec. 3.6). However, the total 
decrease in energy of matter at ~ bpo due to pion condensation, estimated by 
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the two groups: ~ 80 by the Japan and ~ 60 MeV/nucleon by Ihinois, are not 
too different. The recent EOS of the Ilhnois group ||7^ contains this energy gain. 

4.3 Hyperonic Matter 

The possibihty of hyperons contributing to the ground state of dense matter 
has been considered since 1959 |8^. The negatively charged S~ and A~ will 
occur when their chemical potential becomes less than + /ig in matter. With 
A18+(5f +UIX* (A18+(5f ) interactions, the fin + /^e exceeds the rest mass (1197 
MeV) of S~ hyperon at density of ~ 2.2 (2.6) pQ. If the S~-nucleon interactions 
are negligible, they will occur in matter via the weak interaction e~ + n — > E~ +z^e 
at these densities mainly on account of the large pe ~ 170 MeV. Similarly the 
neutral hyperons such as the A will occur when their chemical potentials become 
less than /i^. For the above interaction models the Hn exceeds the A rest mass 
(1116 MeV) at densities of 2.9 and 3.7 po, while the predicted central densities 
of 1.4 Mq stars are 3.4 and 5.1 pQ. 

Much less data exist on hyperon-nucleon (YN) interactions than on the NN, 
and therefore their models are less constrained. LOBHF calculations [|109|| using 



the older Nijmegen soft core YN | 110 and either Paris or Argonne V14 NN inter- 
actions show that the thresholds for and A to appear are not much moved 
by YN interactions. For example, with the Paris NN interaction the threshold 
densities with (without) YN interactions are 3.0 (2.7) and 3.6 (3.7) po for S~ 



and A respectively. The Nijmegen group |111| , |112|] has recently constructed bo- 
son exchange interaction models based on YN and NN data base using SU(3) 
symmetry. LOBHF calculations [|113|| using these models for NN as well as YN 
interactions give 2.2 po for the threshold density of S~, however for A it is pushed 
beyond 7 po- 

The above LOBHF calculations use only two-body interactions without boost 
corrections and TNI, and fail to explain the saturation density of SNM as men- 



tioned in sec. 3.4. Both the BHF groups [ 109| , 113| find that the threshold densities 



are lowered after including TNI effects. For example, the Catania group predicts 
them to be ~ 2.1 and 2.4 po with models including TNI and also with Dirac- 
Brueckner calculations. However, T,~NN and ANN three body forces should 
also be included along with the TNI for consistency. The binding energies of 
A-hypernuclei |114] suggest that the ANN interaction is as strong as the TNI, 



while there is no data on the Ti~NN. 

Ignoring three-body forces and boost interactions in both nucleon and hyperon 
matter, the LOBHF calculations with the unified NN, YN and YY interaction 



model [113] indicate that admixtures of S~ lower the NM energy by only ~ 25 
MeV per nucleon at 5 pQ. At this density the A18-|-(5f-|-UIX* model gives an 
energy of ~ 200 for /3-stable NM. If the energy gain due to T,^ , is not much 
changed by three-body interactions, we should expect ~ 10 % effects on the EOS 
of neutron star matter at p > 2pQ due to hyperons; however, if present they would 
lower the pe and increase the proton fraction Xp significantly. 

4-4 Quark Matter 

When matter is compressed to densities so high that the quark cores of nucleons 
overlap substantially, one expects the nucleons to merge and undergo a phase 
transition to quark matter (QM). The EOS of both HM and QM are necessary 
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to calculate effects of this transition in neutron stars. The interface properties 
are also needed to study the important mixed phase regions. 

At present, lattice QCD can only treat the case of zero baryon chemical po- 
tential and is therefore not useful for neutron stars. Lacking a full theory, the 
simple Bag model is used to estimate the QM EOS. In this model the QM energy 
has a volume term denoted by the bag constant B; it represents the difference 
in the energies of the vacua occupied by hadron and QM, and is responsible for 
the confinement of quarks within nucleons in nuclei. The term dominating at 
high densities is the energy of noninteracting u, d and s-quarks; it is calculated 
neglecting the mass of u, and d-quarks, and using a typical value of ~ 150 MeV 
for that of the s-quark. Since the quark Fermi energies are much larger than 
their masses, the QM properties are not too sensitive to the chosen mass for the 
s-quark. There is no one-gluon exchange interaction energy between quarks of 
different flavor, while that between quarks of of the same flavor i is given by 
{2a/3TT)Ei per quark of flavor i |115f| . Here Ei is the average kinetic energy per 
quark, and a is the strong interaction coupling constant, assumed to have a value 
of ~ 0.5. All higher order gluon-exchange interactions are neglected; their contri- 
bution is presumably subsumed in the bag constant whose value is poorly known. 
Two representative values are B = 122 |ll|] and 200 MeV fm-^ [|T7|. 



The equilibrium conditions for uniform QM containing u, d, s-quarks, electrons 
and muons are: 

fJ-u + fJ'e = fJ-d = /^s, fJ'^l = ^i■e■ (20) 

The energy densities of charge neutral QM and NM are compared in Fig. 12. In 
the interesting region of p ~ 1 fm~'^ the total energy density of quark matter is 
about 1200 MeV fm"^^ ^Nh\c\\ only 122 or 200 MeV fm"^ comes from the bag. 
Assuming the A18+(^t;+UIX* EOS for NM the first order phase transition to QM 
is indicated at ~ 1 baryon/fm'^. 

It has been recently suggested that at high densities QM may have color su- 
perconductivity resulting from non-perturbative attraction between quarks. In 
QM with only u and d-quarks this invariably leads to the possibility of a diquark 



condensate which breaks global color invariance [118|. The associated color gap 



is estimated to be of the order of 100 MeV. When the s-quarks are included there 



are many possible phases |11£]. 



5 Mixtures of Phases in Dense Matter 

The phases of matter considered in the past sections are uniform and locally 
charge neutral, whereas bulk matter needs only to be charge neutral on average. 
For example, iron metal has positively charged regions occupied by iron nuclei, 
and the space in between is negatively charged by electrons. Generally the ground 
state of matter can have a mixture of regions occupied by different phases with 
the constrain of overall charge neutrality. 

Matter in the outer crust of neutron stars, at p < 0.002p0) is believed to be 
like terrestrial matter made up of neutron rich nuclei in electron gas. The inner 
crust, on the other hand, is believed to have a mixture of regions with positively 
charged NM composed of neutrons and protons and PNM with only neutrons 
[p^ j, both immersed in a nearly uniform electron gas. These mixtures occur in 
the density range of ~ 0.002 to 0.6 /?o, beyond which uniform NM is believed to 



be the ground state. More recently mixed phases of NM and QM [12C, 121], and 
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condensates of pions, kaons |122| ], or hyperons in NM have also been considered. 

The Coulomb energy of a single phase uniform matter, due to fluctuations in the 
electron and hadron or quark densities, is negligible; however that of matter with 
mixed phases is not. For example, matter at p = 0.3pQ has drops of NM in PNM. 
Their size is determined by a competition between Coulomb and surface energies; 
large drops have too much Coulomb energy and small drops to much surface 
energy per nucleon. It is necessary to know the energy of the interface to predict 
the nature of the mixed phase region quantitatively. Those for the interface 
between NM and PNM have been calculated from energy-density functionals of 



NM |1C]; while those for interfaces of NM and QM are not well estimated. In this 
section we discuss QM and NM mixed phases, to avoid duplicating the review 
by Pethick and Ravenhall |lO| ], and because they can influence the mass limit of 
neutron stars. 

5. 1 Equilibrium Conditions for Coexistence of QM and NM 

Neglecting surface and Coulomb effects, the equilibrium conditions for the coex- 
istence of of QM and NM at zero temperature are that they have equal pressures, 
and it costs no energy to convert a neutron or a proton in NM into quarks in 
QM. The last condition amounts to 

Pn = 'ipd + Pu and pp = fia + 2/i„ . (21) 

The electron density, and hence the pe, is assumed to be the same in QM and 
NM, therefore the /3-equilibrium conditions are: 

Pn = Pp + Pe and Pd = fJ-s = Pu + Pe ■ (22) 
And finally, the charge neutrality condition is given by: 

fPQM + (1 - f)PNM = eiPe + P^) . (23) 

Here p'qmnm ^"^^ charge densities of QM and NM, /?e,^ are the electron 
and muon densities, and / is the fraction of space filled by QM. A graphical 
representation of these equilibrium conditions is given in ref. |73]. 



As was first pointed out by Glendenning |12C], the QM and NM mixed phases 
can span a wide density region. For example, with the Al8-t-(5f +UIX* NM EOS, 
the mixed phases begin at 4.6 (3.4) and end at 11.3 (9.1) pQ for B = 200 (122) 
MeV fm~^ (Fig. 12). At lower densities we have uniform NM and QM at higher. 
The QM fractional volume increases almost linearly from to 1 in this transition 
range over which /ig — > 0. The uniform QM has almost equal number of u, d, 
and s-quarks, and is essentially charge neutral without electrons. An important 
consequence of /3-equilibrium is that the QM is negatively charged at the beginning 
of the transition, where pe is large. By immersing negatively charged drops of QM 
in the positively charged NM we can remove some of the high energy electrons 
and increase the proton fraction in NM. Both of these lower the energy density 
of matter. 

5.2 Structure of Mixed Phase Matter 

Matter with mixed phases has additional structure due to interfaces dividing the 
regions occupied by the two phases. The surface and Coulomb energies associated 
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with these interfaces, neglected in above section, raise the energy of the mixed 
phase matter as well as determine the topology and structural length scales. The 
surface energy can presumably be estimated from the surface tension, a. For NM 
a ~ IMeV/fm^ whereas for QM it poorly known; typical values are in the range 
10-100 MeV/fm^ [T2l| ). 



Denoting the dimensionality of the structures by D (D = 3 for droplets and 
bubbles, D = 2 for rods and D = 1 for plates) the surface energies are generally 

m 

£s = Da—R" , (24) 

and, for D = 3 the Coulomb energy, £c = {3/5)Z^e^/R, where R the size of the 
structure and 

Y 

is the excess charge of the droplet compared with the surrounding medium. Gen- 
eral equations for £c{R,D) are given in ref. |]lO|. Minimizing the energy density 
with respect to R we obtain the usual result that £s = 2£c at equilibrium. Min- 
imizing with respect to the continuous dimensionality as well thus determines 
both R and D. For droplets {D = 3) the equilibrium radius is found to be: 



Ze = ip'Q\,-p'}!}M)— , (25) 
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A droplet of symmetric NM in vacuum has a surface tension a = 1 MeV-fm~^ for 
which dl^) gives i? ~ 4 fm, which agrees with the fact that nuclei like ^^Fe are 
the most stable form of matter at low density. For QM droplets both the surface 
tension and charge densities are larger but the estimates of R are similar. 

At the beginning of the mixed phase region we expect that spherical droplets 
of QM, with ~ 5 fm, will form a BCC lattice in uniform NM. They would have 
baryon number of the order of few hundred, and a negative charge of similar 
magnitude. As density increases and / approaches 0.5, the drops would merge 
and form rods, which merge further on to form sheats. When / > 0.5 the NM 
sheats break up into rods, and then into drops and eventually disappear when 
/ = 1. This scenario is similar to that in the inner crust; at low densities there 
are drops of NM occupying a small fraction of space. By p ~ O.Qpo NM occupies 
all space via a similar set of mixed phases. 

An other effect of the Coulomb and surface energies is that they decrease the 
density range covered by the mixed phase region. In particular, the lower density 
edge of this interesting region may be pushed up by almost po if a is in the 10 to 
50 MeV fm~'^ range [121]. The energy density of the mixed phase matter is also 



raised by a few MeV fm~^ in this case. Finally, if a were to be large (^70 MeV 
fm~^) the mixed phases may not be energetically favorable, and there will be a 
simple first order phase transition from NM to QM with a density discontinuity. 
One should bear in mind that even if the droplet phase were favored energetically, 
it may not be realized in practice if the time required to nucleate QM drops is 
too long compared to pulsar ages. 



27 



Table 2: Properties of maximum mass (Mm) and 1.4 Mq neutron stars in Mq, 
Pq and km. 



Interactions 



Calc. Ref. p,(M„) ii(M„) p,(1.4) i?(1.4) 



A18 

Al8+(5t; 
Al8+<5t;+UIX* 
A18+UIX 
Paris+TNI 
Bonn A 



Var. 
Var. 
Var. 
Var. 
BHF 
DBHF 



[72] 
[73] 
[73] 

[13] 
[35] 

m 



1.67 
1.80 
2.20 
2.38 
1.94 
2.10 



11.1 
9.4 
7.2 
6.0 
8.3 
6.7 



8.1 

8.7 

10.1 

10.8 

9.5 

10.6 



7.0 
5.1 
3.4 
2.9 
4.0 
3.1 



8.2 
10.1 
11.5 
12.1 
11.1 
11.7 



6 Neutron Star Observations and Predictions 

The gross structure of neutron stars has been predicted using very many EOS, 
phenomenological as well as based on realistic models of nuclear forces [^, 51 1 . Of 



these we consider only those based on realistic models primarily because one can 
always find phenomenological energy density functionals or Lagrangians which 
reproduce their EOS. 

Typical results for nonrotating stars with maximum mass, and with M = 
1.4M0, obtained by recent calculations, are listed in Table 2. The results for A18 
without boost correction 5v are listed primarily for reference. This correction 
is unambiguous [^], and must be added to obtain reliable results. Those for 
A18+5v are also to be taken less seriously, because it gives too large value for 
Pq. The TNI used with the Paris NN interaction ] |85| ] is of the Urbana form 
with parameters determined by reproducing the empirical SNM properties (see 
sec. 3.4). In A18+UIX and Paris+TNI models the 5v is not considered explicitly; 
it is approximately subsumed in the TNI fitted to data. Of the three Bonn 
models, Bonn A comes closest to reproduce the empirical properties of SNM j^] 
with Dirac-Brueckner (DBHF) method. These calculations include the 6v as well 
as many-body forces generated via Z-diagrams. 

The Al8+6v+mX*, A18+UIX, Paris+TNI and Bonn-A (DBHF) models come 
close to reproducing the empirical po; the later two fit the SNM binding energy; 
while the former models fit binding energies of light nuclei via exact calculations, 
since the energy of SNM can not yet be calculated reliably. Nevertheless these four 
^^realistic" models of NM give rather similar results which are not too different 



from those of the 1988 calculations of Wiringa, Fiks and Fabrocini j91] with the 
older Urbana- Argonne interactions now replaced with A18 and UIX. 

The effect of the possible appearance of QM drops in high density matter has 
been studied with the Al8+(5z;+UIX* model. The Mm is reduced to 2.02 and 
1.91 Mq for bag-constant values B = 200 and 122 MeV Fm~^ respectively, while 
the predictions for 1.4 Mq stars remain unchanged. Presence of either kaons or 
hyperons in dense matter is unlikely to have much of an effect on the 1.4 Mq 
stars due to their low central density, while that on the mass limit is difficult 
to estimate quantitatively. For example, if kaons were to condense in matter at 
p = 5po and limit pc to < 5 po, the Mm of the four realistic models will drop to ~ 
2.0, 2.3, 1.7 and 2.0 Mq respectively; while if hyperons were to lower the energy 
of matter at p = 5po by 25 MeV per baryon, the Mm would be reduced by ~ 0.2 

Mq. 

The mass radius relation obtained with models based on the A18 interaction 
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are shown in Fig. 13. Results of A18 and A18+6v models are given primerily 
for comparison. As expected the harder EOS give larger and predict larger 
radii. The differences between the radii predicted by the realistic models is only 
~ 10 %. 

6.1 The Mass Limit 

The observed mass of Hulse-Taylor pulsar B1913+16 of 1.4411 ± 0.00035 
shows that Mm > 1.44M0. All the radio pulsars in known neutron star, and 
neutron star- white dwarf binaries have masses with lower limits less than 1.44M0. 
The X-ray pulsar Vela X-1, which orbits a supergiant, however is consistently 
estimated to have a larger mass of ~ 1.9. The motion of this star is perturbed 
from being pure Keplerian, presumably by tidal forces exerted by the neutron 
star, and its present mass estimate, 1.87to:i7 [0], indicates that Mm > I.TMq at 
95 % confidance level. Finally, if the QPO's indeed originate from the innermost 



stable orbit |18, 19], then Mm > 2Mq. These mass limits are compatible with 
predictions of realistic NM models. 

On the other hand there is no evidence that SN 1987A produced a neutron 
star. Its observed luminosity is now well below the 10'^^ ergs/s Eddington limit. 



suggesting that no neutron star was produced in this supernovae |123]. If we 
assume that the total mass, Mtoi, of the collapsed core plus the matter that fell 
back on to the core after the explosion, went into a black hole, then neutron 
star Mm must be less than ~ Q.S^Mtoi- The factor 0.9 takes into account the 



10 % gravitational binding energy of the neutron star. Bethe and Brown |23] 
estimate Mtoi ~ 1.73M0 using supernovae calculations by Wilson and Mayle, 
and conclude that Mm < 1.56Mq. Uncertainties in these arguments have been 



discussed by Zampieri et. al. [124]. If the conclusion is found to be valid, then 
there must be other explanations for the Vela X-1 observations and the QPO, 
and the NM prediction for the Mm is too large. 

6.2 Temperatures, Cooling and Radii 

Neutron stars are born with interior temperatures of the order 10^^ K, but cool 
rapidly via neutrino emission to temperatures of the order 10^*^ K within minutes 
and ^10^ K in 10^ yr. Spectra observed in X-ray or UV bands for nearby pulsars 
have in some cases black-body components from which surface temperatures of 
order T ~ 10^ K are extracted for pulsars of age 10^ — 10^ years. It is, however, 
unclear how much of the observed radiation is due to pulsar phenomena, to a 
synchrotron-emitting nebula or to the neutron star itself. In other cases upper 
limits have been set from the absence of X-rays. The surface temperatures are 
compatible with predictions from standard modified URCA cooling processes 
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n + n^n+p + e + Ve-, n+p + e -^n + n + Ve- (27) 

Faster cooling processes as direct URCA or due to quark matter, kaon or hyperon 
condensates generally lead to considerably lower temperatures ] |126[ ]. To be con- 
sistent with observed surface temperatures the exotic coolant can only exist in a 
minor portion of the neutron star or it is superfluid whereby cooling is suppressed 
by factors of oc exp(— A/T), where A is the pairing gap. 

The Hubble Space Telescope (HST) has observed one thermally radiating neu- 



tron star RX J185635-3754 with surface temperature T ~ 6 x 10^ K~ 50 eV |20] 
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Its distance is less than 120 pc from Earth and should soon be determined more 
accurately by HST parallax measurements. Circumstancial evidence indicate a 
distance of ~80 pc [51| which leads to a black-body radius of 12 — 13 km from 



its luminosity and temperature. Such radii would agree well with predictions of 
realistic NM EOS (Fig. 13) for M ~ 1 - 2Mq. 

6.3 Glitches and Superfluidity 

Sudden spin jumps, called glitches, superimposed upon otherwise gradual spin 



down have been observed in most of the younger isolated pulsars |127|. Since 
their discovery the Crab and Vela pulsars have each produced about a dozen 
glitches with period changes AP/P of the order of 10~^ and 10~^ respectively. 
In post-glitch relaxation most of the period increase AP decays. 



Many mechanisms have been proposed to explain the glitches [125]. The most 



plausible of these attributes glitches to the angular momentum stored in the 



rotating superfluid neutrons in the inner crust |12g| , 10] . The magnetic torque 
slows down the crust and most of the star except for these superfluid neutrons. 
Their angular momentum is stored in vortices pinned to nuclei in the inner core, 
until an instability occurs that leads to vortex depinning and sudden angular 
momentum transfer to the crust, leading to the glitch. At subnuclear densities 



in the crust, ^Sq pairing between neutrons leads to gaps of order ~1 MeV |125]. 



In NM at p > pq this pairing gap vanishes, but •^P2 pairing of neutrons and S{ 



pairing of protons may occur |128| , |84|. 



The size of the glitches sets a lower limit on the moment of inertia of the 
superfluid in the inner crust which in turn sets a lower limit on the neutron star 
radius for a given mass |127]. Assuming that the mass of Vela pulsar is IAMq, 



a conservative limit on its radius is > 9 km; it is compatible with predictions 
of most EOS. 



7 Conclusions 

Since the discovery of pulsars a significant effort has been devoted to accurately 
calculate properties of dense NM from realistic models of nuclear forces. Exact 
calculations of NM are still out or reach, however the new AFDMC methods 
(sect. 3.7) may eventually succeed. The present variational upper bounds seem 
to be above the true energies by ~ 12 %. Such an error does not have serious 
consequences on the predicted properties of neutron stars. For example, an EOS 
obtained by reducing the variational energies, without rest mass terms, by 12 % 
reduces the maximum mass of AlS+Sv+VTK* model by 2.3 % to 2.14 Mq, and 
the radius of 1.4 Mq star by 2.9 % to 11.2 km. Larger uncertainties stem from 
the fact that the double vr'^ and vr" condensation scenario illustrated in Fig. 11 
has not yet been calculated with realistic interactions, though it appears unlikely 
that it will influence the NM EOS by much more than 10 %. 

Local models of two-nucleon interaction seem to be now converging. The pre- 
dictions based on the 1988 calculations with Argonne 14 interaction are not too 
different from those of the 1998 calculations with the more accurate A18. It also 
seems likely that the local models give a fairly accurate description of two-nucleon 
interaction. A concern is that the present models of TNI are based on fits to a 
rather limited set of data, and are not as precise as the NN-interaction models. 
However, addition of the UIX* TNI to the A18+6v increases the maximum mass 
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by ~ 20 % and -R(1.4) by 13 % (Table 2). These changes may be important but 
they are not very large. 

The present models of kaon-nucleon and hyperon-nucleon interactions are based 
on very limited data, and we have none on K~ NN and Ti~ NN three-body forces. 
These could have significant effect on the threshold densities for kaons and hy- 
perons to appear in dense matter. Hopefully advances in QCD and quark-models 
will provide a more rigorous framework to describe these interactions, and calcu- 
late properties of quark matter. The bag model estimates of QM EOS may have 
significant corrections at densities of interest in neutron stars. 

^Prom present observations there seem to be three possible scenarios for the 
limiting mass of neutron stars. If QPO's are indeed due to accretion from the 
innermost stable orbit, then the NM predictions of ~ 2.2M0 are reasonable, 
and strange baryons and quark drops do not soften the EOS of matter at /) < 7po 
significantly. If the Vela X-1 mass measurement is correct, but QPO's have some 
other origin, then Mm could be ~ 1.8 M©, indicating some softening of the 
NM EOS. However, if the present interpretation of QPO's and Vela X-1 mass 
measurements are both faulty, and Mm is as small as 1.56 Mm as estimated from 
the absence of a neutron star in SN 1987A, then a significant softening of the NM 
EOS by phase transitions is indicated. Further observations will hopefully clear 
this situation. 

Phase transitions such as NM to QM, can soften the EOS significantly. Fortu- 
nately these can have a measurable effect on the spin down of a rapidly rotating 



star in favorable cases, as has been recently pointed out |130, 131 1. Consider the 
case of a rapidly rotating star whose central density is close to a first order phase 
transition. As the star slows and the central pressure increases due to decrease 
of the centrifugal force, the core matter will change its phase and become more 
dense at a critical angular velocity Oc- This decreases the moment of inertia, 
which assumes the characteristic form: 

I = lJl + Cif]2 _ C2(f)2 _ f]2)3/2 ^ \ _ ^28) 



for < f^c- Here, ci and C2 are small parameters proportional to the density 
difference between the two phases, and C2 = for > ric- 

In order to make contact with observation, the temporal behavior of angular 
velocities must be considered. The pulsars slow down at a rate given by the loss 
of rotational energy, believed to be given by: d{^IO,'^)/dt oc — fi^^^, where n = 3 
for dipole radiation, Eq. (|^) and n = 5 for gravitational radiation. With the 
moment of inertia given by Eq. (^) the angular velocity can be calculated. The 
corresponding braking index, n{Q) = Ofi/O^, depends on the second derivative 



of the moment of inertia, I" = dl/d fl. Using Eq. (28) we obtain: 



n{n) ^ n - cin^ + C2-^^=^ . (29) 

which exhibits a characteristic (l^c — ri)~^/^ singularity as 0, approaches from 
below. Observations of the braking index of a rapidly rotating, new born pulsar 
would be very interesting. 

All realistic NM EOS predict that the radius of neutron stars with a mass of 
1 to 1.5 Mq is ^ 11 to 12 km. Future high resolution Chandra and XMM space 
observatories will hopefully be able to measure black-body spectra and detect 
gravitationally redshifted spectral lines from several stars. Such observations will 
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help determine masses, radii and temperatures uniquely if the distance of the 
star is known. It is important to know the radius of a 1.4 Mq star, because that 
would test the EOS in the p^3pQ region in which large modifications of NM EOS 
are not expected on the basis of our present, naive estimates of kaon-nucleon and 
S~-nucleon interactions. 
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Figure 1: The energy per nucleon in uniform neutron matter at low densities, 
comparison of various effective interactions: Reinhard et al. 
~ ~ SkM* Nijmegen NL3 and NLZ ] 
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Figure 2: Deuteron radial wave-functions, u(r) (upper curves) and ^Di w(r) 
(lower curves) predicted by modern NN interactions: CD-Bonn (solid), Nijm-I 
(dashed), Nijm-II (dash-dotted), Reid 93 (dotted) and A18 (long-dashed). 
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Figure 4: Static part of the NN potential in the deuteron in spin projection 

M = state. The upper curves show the potential along the Z-axis (9 = 0), 
while the lower curves show it in the X — Y plane. The '^0 denote the central 
and tensor components of the NN interaction in the deuteron. 
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Figure 5: Charge densities of two protons located one fm apart at Z = ±0.5, 
obtained by inverting tlic dipolc approximation to proton charge form factor. The 
sharp peaks at Z = ±0.5 arc unphysical, they will be rounded off by relativistic 
corrections and improved data on proton form factor. 
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Figure 6: The observed energies of all bound and quasi-bound states of up to 
eight nucleons are compared with the preliminary results of GFMC calculations 
with < 2% errors, using A18 model of Vij and Illinois model 2R of Vijk 
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Figure 8: Comparison of the energies of PNM and SNM obtained for the A18 
model with the variational chain summation (VCS) method [^3| and LOBHF 
(LOB) method |83|. The true results for this interaction are expected to be few 
MeV below the VCS. 
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Figure 9: The PNM and SNM energies for the A18+5v + UIX* model, and the 
fits to them using effective interactions. The full lines represent the stable phases, 
and the dotted lines are their extrapolations. From [ |73| . 
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Figure 10: Kaon energy as function of NM density calculated with the various 
approximations discussed in the text is compared with the electron chemical 
potential //e calculated from A18+(5?;+UIX* model in Ref. [ffsf . 




Figure 11: Schematic drawing of the spin arrangement of neutrons and protons 
in a phase with vr'^ condensation. The nucleons are expected to reside mostly in 
nodal plains of the tt^ field, where the field gradient is largest. The charged vr" 
may condense with momenta perpendicular to that of the vr*^ field. 
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Figure 12: The energy densities of NM (full lines), QM (dashed lines) and 
mixtures (dotted lines) from [73|. The rest mass contribution m„p to the en- 
ergy density of NM is subtracted from the results of all the models for easier 
comparison. 




Figure 13: Neutron star gravitational mass, M(i?), in solar masses vs. radius in 
kilometers for the four models described in the text. Full curves are for /3-stable 
matter and dotted ones for pure neutron matter. The dashed curve, FPS, is from 
H]. Fig. from @. 



